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ABSTRACT 



Context. Several studies have investigated the fractal and multifractal nature of magnetic features in the solar photosphere and 
its variation with the solar magnetic activity cycle. 

Aims. Here we extend those studies by examining the fractal geometry of bright magnetic features at higher atmospheric levels, 
specifically in the solar chromosphere. We analyze structures identified in CallK images obtained with the Precision Solar 
Photometric Telescopes (PSPTs) at Osservatorio Astronomico di Roma (OAR) and Mauna Loa Solar Observatory (MLSO). 
Methods. Fractal dimension estimates depend on the estimator employed, the quality of the images, and the structure identifi- 
cation techniques used. We examine both real and simulated data and employ two different perimeter-area estimators in order 
to understand the sensitivity of the deduced fractal properties to pixelization and image quality. 

Results. The fractal dimension of bright 'magnetic' features in CallK images ranges between values of 1.2 and 1.7 for small and 
large structures respectively. This size dependency largely reflects the importance of image pixelization in the measurement of 
small objects. The fractal dimension of chromospheric features does not show any clear systematic variation with time over the 
period examined, the descending phase of solar cycle 23. 

Conclusions. These conclusions, and the analysis of both real and synthetic images on which they are based, are important in 
the interpretation of previously reported results. 

Key words. Sun: activity - Sun: chromosphere - Sun: faculae, plages 



1. Introduction 

The spatial distribution of the solar magnetic fields is very 
complex, depending on both the varying level of mag- 
netic activity over the solar cycle and the height of ob- 
servation in the solar atmosphere. This complexity likely 
results from both the dynamo process itself, which may 
occur on many different spatial scales, and the interac- 
tion of the field with convective motions as it emerges 
through the Sun's outer layers. The signatures of these 
processes have been investigated in previous works by frac- 
tal analyses of solar active regions, but the quantitative 
results obtained differ widely depending o n the type of 
data an d analysis techniqu es employed fe.g. lJanfien et all 
[2003: McAteer et al.Ll2005l) . Moreover, the lack of a unique 
definition of the fractal dimension itself often makes com- 
parison of results difficult. 
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Among recent studies, one focused on the possible re- 
lationship betwe en magne t ic fea ture complexity and so- 
lar cycle phase l|Meunieil l2004j) . This study looked at 
a large sample of active region magnetograms acquired 
with the Michelson Doppler Interferometer (MDI) aboard 
the Solar and Heliospheric Observatory (SOHO) between 
April 1996 and June 2002. The measured fractal d imension 
increa sed with structure size ( i n agre ement with iMeunieJ 
l|l999|l and lNesme-Ribes et all ()l996|) '). showing a peculiar 
change in behavior near structures of area 550-800 Mm^. 
A similar dependen ce on structure size was also found by 
Ijanfien et alJ l)2003jl in both high resolution photospheric 
magnetograms and r nagnetohydrody namic (MHD) simu- 
lations. Additionally, iMe^mieil (|2Q^) investigated the re- 
lationship between the geometry of facular structures of 
different spatial scales and magnetic field intensity, flare 
activity, and solar cycle phase. She found that, while a 
region complexity generally increases with magnetic field 
intensity, there is no clear correlation with flare activity. 



2 



S. Criscuoli et al.: On the reliability of the fractal dimension measure of solar magnetic features 



Variations of fractal dimension with solar cycle were also 
reported, but their amplitudes and sign largely depended 
on object size and the associated magnetic field. 

In order to investigate the complexity of magnetic fea- 
tures using observations representative of chromospheric 
heights, we have analyzed the fractal dimension of bright 
features identified in full-disk CallK images acquired by 
the Precision Solar Photometric Telescopes (PSPTs) at 
Osservatorio Astronomico di Roma (OAR) and Mauna 
Loa Solar Observatory (MLSO). The data analyzed span 
the past 6 years and thus allow investigations of variation 
with the solar cycle. 

Several factors can influence fractal dimension estima- 
tion. Both image resolution and projection effect corre- 
lation, mass function, and perimeter-are a estimators in 
studi e s of interstellar mol e cular clouds ijSanchez et al.l 
I2nn5t IVogelaar fc Wakkeil llQQ^ - and resolution and 
thresholding effects are important in fractal dimension 
estimation of flow patte r ns in fi eld soil by box-c o unting 
methods ()Baveve et al.l Il998|) . Lawrence et all ()l996j) 
studied similar effects in multifractal and fractal measures 
(box-counting, cluster dimension, threshold set) of solar 
magnetic active regions. The effec t of struc t ure selectio n 
technique was also investigated hv lMeimieT<lll999H2nn4) . 

To link our findings dir ectly to those o f several recent 
solar studies (for example iMeunierl l2004t Ijanfien et al.l 
we have employed the perimeter-area relationship 
to define the fractal dimension of the identified features. 
To interpret our results, we have investigate the sensitivity 
of the deduced fractal dimension to the pixelization and 
resolution of the image and to the perimeter measure al- 
gorithm employed. In particular, we have determined how 
these factors influence the geometric properties deduced 
for objects of different sizes, and more generally have ad- 
dressed the question of whether the perimeter-area rela- 
tion is suitable to the study of the fractal and multifractal 
nature of solar magnetic features as a function of the solar 
cycle. 

This paper is organized as follows. In the next sec- 
tion we describe the observations, data processing tech- 
niques, and geometrical measures employed. In section §3 
we present the results obtained and in §4 compare them 
to those of previous efforts. In §5 we investigate the sen- 
sitivity of the deduced fractal dimension to the image 
resolution and the measurement techniques employed, by 
examining synthetic structures whose fractal properties 
are theoretically known: non-fractal objects, von Koch 
snowflakes, and those produced by fractional Brownian 
motion (fBm) . The effect of atmospheric seeing is also in- 
vestigated using PSPT images taken under variable seeing 
conditions. In section §6 we discuss our and previous re- 
sults, in light of the conclusions drawn in §5. Finally, in §7 
we summarize our work and conclude with a more general 
discussion on the validity of the methods. 



2. Observations, processing and definitions 

2.1. PSPT data 

The bulk of the data we analyzed is from the archive of 
daily full-disk observations carried out with the PSPT at 
OAR. This was supplemented with data from the PSPT 
at MLSO for consistency and resolution tests (§3.1 and 
§5.2). Details abo ut the data an d the im age pre-processing 
can be found in lErmolli et alJ l|200(t l. In brief, the im- 
ages were taken with "twin" telescopes at the two sites, 
through interference filters centered at three wavelength 
bands (CallK line center 393. 4nm, fwhm 0.27nm, blue 
continuum 409. 4nm, fwhm 0.27nm, and red continuum 
607. Inm, fwhm 0.46nm), with a 2048 x 2048 16 bit/pixel 
CCD camera, yielding a spatial scale of ^ 1" per pixel. 
Images from OAR were binned to half resolution, yield- 
ing a final spatial scale of ^ 2" per pixel. Images from 
the two telescopes were independently dark and flat-field 
corrected and had the mean center-to-limb variation re- 
moved. The images of any one wavelength triplet were 
resized and aligned to allow pixel by pixel comparison be- 
tween filters. 

For this study we selected OAR daily image triplets, 
obtained on 238 different observing days during the sum- 
mers (July to September) of 2000 through 2005. We chose 
images acquired during the summer months because these 
are generally of higher quality. In order to compare results 
obtained with the two instruments, we also selected 44 
triplets (the best in the CallK band of the day, according 
to quality criteria described in §2.2) from the MLSO and 
OAR archive taken during the summer of 2005. For that 
comparison, MLSO images were rescaled to match OAR 
spatial scale images (~ 2" per pixel). 

Finally we were able to quantify the effects of atmo- 
spheric seeing by using MLSO images acquired at 10 min- 
utes intervals throughout the day, weather permitting. 
According to the quality criteria explained in next para- 
graph, we selected 27 pairs of high and low quality triplets 
(one pair per day) from the period February to October 
2005. For this analyses, the full resolution 1" per pixel) 
MLSO data were employed. 

2.2. Data quality 

The geometric properties of solar features extracted from 
the images are likely sensitive to the spatial resolution of 
the image being analyzed. This in turn depends on atmo- 
spheric and instrumental operation conditions during the 
observation. To estimate the inherent quality of any given 
data image, we measured (in pixel units) the width of a 
Gaussian fit to the limb profile observed in CallK images. 
Small values of the solar limb width indicate lower instru- 
mental or atmospheric smearing and thus better quality 
images. The limb width distributions of our datasets are 
asymmetrically shaped with a long tail toward higher val- 
ues, so that the mean is not the most probable value. The 
mean limb profile width of the OAR CallK (binned to 
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groups, the best and the worst images were selected, so 
that for each observing day, two different quality images 
were retained. This restricted our analysis of seeing effects 
to 27 observing days. The mean limb width for the two 
groups were 3.5±0.2 and 4.9±0.3 for high and low quality 
images respectively. 



Fig. 1. Central disk detail of a CallK image from OAR 
PSPT archive and corresponding mask obtained with 
the second identification method described in the text. 
Features with areas lower and larger than 2000Mm^ are 
shown in gray and black colors respectively. 



half resolution) images analyzed is 2.5 ± 4.0 pixels with a 
median value of 1.4 pixels. That of the MLSO CallK im- 
ages from the summer 2005 is 4.1 ±0.8 pixels for the mean 
and 3.9 pixels for the median (measured on full resolution 
images), while that for the OAR images acquired in the 
same period are 4.0 ±6.5 pixels and 2.2 pixels respectively 
(measured on binned half resolution images). Considering 
the different pixel scale of MLSO and OAR images, the 
two 2005 datasets have similar median quality (about 4 
arcsec), but the OAR dataset contains a higher number 
of low quality images, skewing the mean to a much higher 
value. 

The MLSO limb width distribution for year 2005 has 
its maximum at a value of about 3.7 pixels (measure on 
full resolution images). To study fractal dimension depen- 
dence on seeing condition (§5.2), images from each day 
were grouped in two sets: those whose limb widths lie be- 
low O.lcr from the peak, and those whose limb widths lie 
between 0.8cr and 3. Oct above the peak. From each of these 



2.3. Feature identification 

Bright features were identified in the CallK images using 
two methods based on a combination of pixel intensity and 
connectivity. Because of the link betw een CallK bright- 
ness and magnetic flux intensitv (e.g. Skumanic h et al. . 
TigiHaTvev fc Whitel.ri99iEir2003: Ortiz and Rast . 



200^ these features are assumed to represent small mag- 



netic structures, but this assumption plays no role in their 
identification. 

The first ident i ficatio n method is analogous to that 
used by iMeunied l|2004j) . taking into account intensity 
alone and bas ed on fixed thresholding values. While 
iMeunieJ l)2004(l employed thresholding on magnetograms, 
we select pixels whose intensity contrast ^ in the CallK 
images exceeds a given value. Umbral, penumbral and pore 
pixels, were included selecting those pixels whose intensity 
in red continuum images was below a given threshold. 

The second identification method applied takes into 
account both pixel intensity and connectivity. As de- 
scribed in detail by lErmolli et all l|2006(l . we first iden- 
tified regions of the solar disk which include active re- 
gions and their remnants. These pre-selected regions are 
made up of those pixels which are brighter than a fixed 
contrast value on a spatially filtered image (box average 
side — ) . A second contrast threshold value, deter- 
mined by the procedure explained in iNesme-Ribes et alJ 
l|l996() . is then used to single out from the pre-selected 
regions all those features identified for our study. Pixels 
representative of sunspots and pores, previously identified 
from red continuum images, were excluded from the iden- 
tified features. 

For the subsequent fractal analysis, the two identifica- 
tion methods were employed to produced two independent 
binary masks from each image triplet processed. In each of 
these pixels satisfying one of the two identification criteria 
above were assigned a value of one, with all other pixels 
set to zero. An example of a CallK PSPT image (showing 
only a central disk region) and the corresponding mask, 
obtained with the second identification method described, 
is given in figQ To reduce distortion due to projection ef- 
fects, the analysis was restricted to structures near disk 
center, /i > 0.8, where /i is the cosine of the heliocentric 
position angle. Additionally, isolated bright points were 
removed from consideration by discarding all structures 
of area less than 10 pixeP. 



Jc = , where 7 is the intensity measured at each pixel 
and Jo is that representative of the quiet sun and obtained by 
a fit to its center to limb intensity variation. 
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2.4. Perimeter and area evaluation 

There are several ways to define and evaluate the 
perimeters and areas of fe atures in a binary image 



I Gonzalez and WoodsL l2002|) . and thus characterize the 
independent structures. The goal is to define, detect, and 
count the pixels which constitute the feature edges. For 
our study we considered three methods, and evaluated the 
errors associated with them. 

In the first method, we defined border pixels by row 
and column, identifying for each the pixel for which the 
binary value changes. The perimeter was then evaluated 
by summing the external sides of the border pixels, so that 
for example an object made up of 1 pixel has an area of 1 
and a perimeter of 4, while one made up of two pixels has 
an area of 2 and a perimeter of 6 or 8 depending on the 
pixels' relative positions. 

In the second metho d, we applied the Roberts opera- 
tor ( Turner et al.L Il998l) to the image in order to identify 
border pixels and defined the perimeter as the sum of the 
all pixels whose value is not zero. Using this method, an 
object of 1 pixel has a perimeter of 4 and an object of 
two pixels always has a perimeter 6, independent of the 
relative positions. 

In the third method, pixels are identified as border 
pixels if they are connected from between 1 and 7 of the 
neighboring 8 contiguous pixels. The perimeter is the sum 
of the selected pixels, so that an object 1 pixel in area has 
a perimeter of 1 and an object of area 2 has a perimeter 
of 2, independent of the pixels' relative positions. 

Only the results obtained using the first method are 
included in the body of this paper. The reasons for this 
are discussed in §5.1. 

2.5. Fractal dimension definition 

Several definitions of the fractal dimension of two- 
dimensional structu res and correspondi ng techniques for 
its estimation exist l|Turner et allll99^ . 

If a structure is self-similar, its perimeter L and area 
A display a power-law relation: 



L oc A"'^ 



(1) 



where d is the fractal dimension. With this definition, 1 < 
d < 2 and d~l for non-fractal structures. 

We estimated d using two methods. In the first, we per- 
formed a simple linear fit to the logarithm of the perime- 
ters and areas measured for structures of different sizes; we 
indicate the fractal dimension so obtained as D. In order 
to investigate the size dependence of the fractal dimen- 
sion estimated in this way, the fit is performed over the 
entire data set or over objects in a specified size range. 
In the second, we adopt ed a method first pro posed by 
iNesme-Ribe s et ahl ()l996|) and later employed bv lMeun'ien 
(199 9, 20QA) , in which perimeter and area values are aver- 
aged over bins in area, each of width A log A = 0.05, and 
the fitting is done on these averages for a series of overlap- 
ping windows of constant width A log A =1.5, producing 
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Fig. 2. Fractal dimension dl versus area of bright features 
identified on calcium images (~ 2" /pixel). Full circles: 
summers 2000-2005 OAR-PSPT. Open triangles: summer 
2005 MLSO-PSPT. Full triangles: summer 2005 OAR- 
PSPT. dl increases fast with object size at area smaller 
than 2000 Mm^. For larger areas, a plateau is observed for 
summer 2005 OAR and MLSO data, and a slow rise on 
the 2000-2005 OAR dataset. 



OAR2000 - 2005 


OAR2005 


MLS'02005 


1.337 ± 0.002 
1.64 ±0.02 


1.307 ± 0.003 
1.53 ±0.07 


1.307 ± 0.004 
1.54 ±0.09 



Table 1. Fractal dimension D estimated for features se- 
lected on OAR-PSPT and MLSO-PSPT CAIIK images. 
First row: fractal dimension D considering the entire range 
of structure areas. Second row: fractal dimension D con- 
sidering only structures larger than 2000 Mm^. 



a measure of d which is a function of A. We indicate the 
fractal dimension estimated in this way as dl. For both 
methods linear fits were performed by a chi square mini- 
mization and the associated error is taken to be the vari- 
ance in the estimate of the slope. 

More details about the significance of the two estima- 
tors and how they relate to each other are given in the 
Appendix. 

3. Results 

3.1. Fractal dimension and feature size 

Figure |21 shows the variation in fractal dimension dl of 
the identified chromospheric features as a function of their 
size, as derived from the OAR and MLSO PSPT CallK 
rebinned data using the second identification method de- 
scribed in §2.3. The results obtained from three data sets 
are shown: the full 2000-2005 OAR summer period, the 
single 2005 summer OAR data, and the single summer 
2005 MLSO data. For all data sets, dl increases with ob- 
ject size, increasing fastest for structures of smallest areas 
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and becoming almost constant at the largest scales. The 
three curves overlap for structures of area less than about 
1000 Mm^, corresponding to about 500 pixel^. For objects 
of size greater than about 1000-1500 Mm^ the 2005 OAR 
and MLSO data both show a plateau in the measured 
fractal dimension. Somewhat surprisingly, this plateau is 
less evident when analyzing structures from the full 2000- 
2005 OAR data set. The fractal dimension deduced over 
this longer period continues to slowly increase even at the 
largest scales. 

Table ^ shows the value of D obtained from the dif- 
ferent data sets when a single fit to the perimeter area 
relation is made over structures of all sizes (top row) or 
only those of area larger than 2000 Mm^ (bottom row), the 
threshold value suggested by the trends observed in figEl 
In agreement with dl estimates, the fractal dimension D 
is reduced by the inclusion of the small and apparently 
less complex regions. In this case, because of the large 
number of objects with small areas, D is biased toward 
a low value and the formal error quoted is quite small 
(as is a measure of the fit) in spite of the fact that a 
single linear fit does not reflect the perimeter area rela- 
tion at all scales (see also Appendix). The large objects, 
who's perimeter area relationship is poorly fit by the single 
slope estimator, are insufficient in number to significantly 
alter the fit value or infiuence the error measure. Figure 
n displays typically structures of areas larger and smaller 
than 2000Mm ^. The smallest objects appear somewhat 
rounder and more regular than the largest ones, but no 
overwhelming difference between the two groups can be 
inferred by visual inspection alone. 



3.2. Temporal variation 

Figure 13 shows the variation of fractal dimension dl with 
both feature size and time for the six year OAR period 
analyzed. For sake of clarity, only the largest error bar 
(belonging to the largest area objects of year 2004) is in- 
dicated on the plot. The fiuctuations in the other values 
shown are generally smaller than the largest differences 
observed among the years. We find that the variation in 
fractal dimension does not show a clear correlation with 
solar cycle over the period analyzed (the descending phase 
of Solar Cycle 23). The values for large structures show 
significant year to year variation, with the maximum and 
the minimum dimensions measured for years 2002 and 
2005 respectively. The reliability of the 2005 values are 
supported by the nearly identical results obtained from 
the independent OAR and MLSO measurements (Table 1, 
fig|2J). The other years show a plateau value at about 1.6, 
although in the area range 2500-7000 Mm^ year 2000 has 
a mean value of about 1.65 and a slight increase with ob- 
ject size up to a maximum value of 1.7 is measured for 
year 2002. Note that if we restrict the analyses to the 
area range 2500-7000 Mm^ then a weak trend with the 
last cycle, showing a double activity peak in 2000 and 
2002, is observed. Figure^ shows the temporal variation 
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Fig. 3. Temporal variation of the fractal dimension dl 
versus area for features identified on OAR-PSPT calcium 
images. The bar on the left represents the largest error 
bar, obtained for the largest areas for year 2004. At area 
smaller than about 1000 Mm^ all the curves overlap, while 
differences (not clearly correlated with solar cycle) are ob- 
served at the largest areas. 
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Fig. 4. Temporal variation of the fractal dimension D ver- 
sus area for selected OAR-PSPT calcium images and for 
different area range. Error bars in the case of fits per- 
formed on the whole dataset (circles) or at smallest ob- 
jects (triangles) are smaller then the symbol size. Results 
obtained for the largest area are in good agreement with 
results obtained by dl estimator (figOJ. 



in D for three different area ranges. For the reasons ex- 
plained in §3.1, error is omitted when fit is performed on 
the entire area range and at the smallest areas. In this 
measure a small trend with solar cycle is observed when 
the entire structure size range is included in the fit. When 
fit is performed on the largest objects, the highest frac- 
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tal dimension is measured for year 2002, and the lowest 
for year 2005, and similar values are measured for the 
other years, in agreement with results obtained for dl. For 
smaller magnetic regions, again the maximum is observed 
in 2002, but a minimum is found for year 2004 data. We 
notice that, while suggestive, the trends at larger objects 
are small and the fractal dimension is constant within the 
measurement uncertainties. Moreover variations of coeffi- 
cients evaluated over the entire area range and at smaller 
areas probably reflect variations in the size distribution o f 
the magnetic regions l)Meunieii 120031 : lErmoUi et alll2006| ) 
rather than a real temporal variation of the fractal dimen- 
sion. 



4. Comparison to previous results 



4.1. Fractal dimension and structure size 



Since the work of iRoudier and Mulled l)l987|) the frac- 
tal geometry of structures found in images of the 
outer layers of the solar atmosphere has been in- 
vestigated by a number of authors. Both magnetic 
features, at moderate to high spatial resolution, and 
non-magnetic features ass ociated with plasra a mo- 
tions have been studied llRoudier and MulleJ. Il98 



Lawrenc e. Ruzmaikin & Cadavi 
Nesme -Ribes ct alj. 



Balke et al 
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199St lMeunie4^ fiflfla IStenflo and Holzreuteri M 
JanBen et all l2003t iMeunieri 120041: iMcAteer et all 12005') 



mien i:^t»)4t IMC Ateer 
In order to ensure a meaningful comparison, we will 
compare the results we have found here only to those 
previously published results for active regions investigated 
using the perimeter-area estimator. 

In agreement with previous results, we find a fractal di- 
mension that increases with feature size, from a minimum 
value of about 1.2, to an approximately constant value of 
1.5 — 1.7 for structure areas larger than ^ 1000 — 2000 
Mm^. This range in dl agrees well with measurements by 
Nesme-Ribes et al ] (|l99ftl . but not with those reported by 



Meunien l|l999t l2004l) , who found generally a higher min- 



imum value (around 1.4). From the analysis of both real 
and simulated data (§5) we know that the minimum value 
measured is somewhat dependent on both the identifica- 
tion method employed and the image resolution. We sug- 
gest that the higher minimum value reported bv iMe^mieil 
l(l999t I2OO4I) may be a consequence of image resolution 
(see §5.3), as the full-disk MDI data she analyzed are un- 
affected by atmospheric degradation. The plateau in dl 
beyond object sizes of 2000 Mm^ also agrees wi th previ- 
ous re sults, but this time more so with thos e of [ Meunierl 
l|2004l) and less so with Nesme- Ribes et al.l l|l99fi() . who 
found the plateau to occur already for structures of size 
> 300 square-pixels (corresponding to about 500 Mm^). 
Tests with both real and simulated data suggest that this 
difference may lie in the area range over which the fit 
for dl at each p oint was made, a value not quoted by 
iNesme-Ribes et al.. (i,1996>) . but taken in our study to be 



A log A — 1.5 in agreement with that used bv iMeunieJ 
(I2OO4I) . 

An increase of fractal dir nension with rnagne tic fea- 
ture size was also observed bv l.Tanfien et all (2Q^. They 
studied fractal dimension of magnetic features analyzing 
high resolution (approximately 0.4") magnetograms 
acquired with the Vacuum Tower Telescope and synthetic 
images obtained through MHD simulations. They found 
D = 1.38 ± 0.07 on synthetic data, and D = 1.21 ± 0.05 
on real data, corrected to D = 1.41 ± 0.05 when taking 
in to account resolution effects. Note that our estimate, 
D = 1.3 (fit on entire area ra nge), lies in between these 
last two values. The plots of Ijanfien et alJ l)2003|) show 
a deviation from these fits for log(A/pixel^)>2.5 (corre- 
sponding to 315 pixel^) for both real and simulated data, 
despite the differing pixel scale of the two data sets (~ 72 
km for real data and ~ 21 km for simulated ones). Fits 
to objects whose areas were larger then this threshold 
gave D = 1.47 and D = 1.9 (values not corrected for res- 
olution effects) for real and simulated objects respectively. 



4.2. Temporal variation 

iMeunieJ l|2004|) performed a time-dependent analysis, 
evaluating the variation in the fractal dimension dl with 
object size for three different periods: minimum, ascend- 
ing and maximum phase of the current solar cycle. A cor- 
relation with solar activity for structures of size ~ 1000 
Mm^ was reported, with the highest fractal dimension be- 
ing measured during the cycle maximum period. Larger 
structures (2000-7000 Mm^) were found to have a higher 
fractal dimension during the ascending phase of the cycle 
than at cycle maximum. Variations were of the order of 
few per cent. The same trends were found for estimates 
of D, but with larger amplitude variations. If we restrict 
our analyses to the area range 2000-7000 Mm^, we instead 
find a little correlation of dl with solar cycle, the highest 
values being measured for years 2000 and 2002, and the 
smallest for 2005. The amplitude of the variatio ns in our 
data i s slightly higher than the one reported bv iMeunied 
1)2004!) . the largest yearly variation measured over the six 
year period being of order 10%. The trend reported for 
structures of moderate size (1000 Mm^) is not observed in 
our analysis. 

5. Discussion of fractal dimension estimation 

Assessment of the fractal dimension of features in digital- 
ized images requires a series of operations: 

— Image segmentation to isolate regions of interest, 

— Edge identification in the resulting bi- level images, 

— Perimeter and area measurement of structures so iden- 
tified, 

— Fractal dimension evaluation using these measures. 

Each of these steps introduces a certain degree of arbitrari- 
ness which influences the result. Moreover, the results are 
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sensitive to intrinsic differences between image sets, un- 
related to the geometric properties of the features they 
capture. In this section we focus on the effects of edge 
identification technique, pixelization, and resolution, by 
analyzing synthetic images of objects whose fractality is 
known: non-fractal objects, the von Koch curve, and ob- 
jects obtained by fractional Brownian motion. Seeing ef- 
fects are also investigated through the analyses of MLSO 
PSPT data. 




-T ' 1 ' 1 ' 1 ' 1 ' r 

O 1000 2000 3000 4000 5000 



Area (pixel ) 



5.1. Perimeter definition and pixelization effects 

To study the influence of the perimeter finding algorithm, 
we examined the empirical dimension of non-fractal ob- 
jects as a function of their size. Three different perimeter 
identification techniques (described in §2.4) were applied 
to three geometric shapes (squares, right triangles, and 
circles). In the absence of error, all three methods should 
yield a value of one since the objects are non-fractal, but 
because of image pixelization, fractal dimensions greater 
or lower than one were measured. 

We found that errors in fractal dimension evaluation 
are functions of the object size for both D and dl. Figure 
|S1 shows the results obtained for circular objects, with the 
top panel showing dl versus object area and the bot- 
tom panel plotting D, evaluated by fitting points of area 
greater then a given threshold, as a function of the thresh- 
old value itself. In both cases, errors are greatest for ob- 
jects of small size but persist to surprisingly large scales. 
Analogous trends were observed for the other shapes ana- 
lyzed. In dl estimations, errors of less than 5% are achiev- 
able for object sizes greater than some hundreds-1000 
pixel^, but for circular objects, which can not be grid 
aligned, the error never drops below 1%, independent of 
the perimeter measure employed. This is true even for ob- 
ject sizes exceeding 5000 pixel^. For any given size object 
D is significantly closer to its expected value of 1 than is 
dl. This is because the evaluation of D in the perimeter- 
area fit is performed over all points above a minimum 
size. This includes the large objects not included at small 
scales in the evaluation of dl. Therefore, the object size 
threshold above which the error in D is below 5% occurs 
at smaller scales than for dl, but still is not usually less 
than some hundred square pixels. 

The origin of these errors lies in the impossibility of 
representing curves or non-grid aligned lines on a rect- 
angular grid. This causes the area and the perimeter to 
scale differently from what is expected for non fractal ob- 
jects. For instance, in the case of a right triangle whose 
two sides are grid aligned, the overestimation of the hy- 
potenuse leads to the overestimation of both perimeter 
and area. It can be shown that, because the relative error 
in the perimeter estimation is not size dependent, while 
the relative error in area estimation decreases with in- 
creasing object size, the estimated fractal dimension is al- 
ways overestimated. 




1000 2000 3000 

Object size threshold (pixel^) 

Fig. 5. Fractal dimension dl (top) and D (bottom) esti- 
mated for circles as function of object size and minimum 
area threshold respectively. Results obtained with external 
sides (crosses) and 8-contiguous point (triangles) perime- 
ter finding algorithm are shown. Error bars are smaller 
then symbol size. Note that the horizontal scales differ for 
the two panels. 
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1.8 


1.516 


16% 


fBm 


1.6 


1.413 


12% 


vonKnoch 


1.26 


1.310 


-3% 



Table 2. Theoretical fractal dimension Dth, measured 
fractal dimension D for the studied objects and the rela- 
tive error e = {Dth — D)/ Dth- Pixelization errors increase 
with increasing structure complexity. 



For our analysis of solar data, we employed only the 
row and column counting method (the first method de- 
scribed in §2.4). This method was chosen because it alone 
produced no error for grid aligned squares, a minimum 
criterion. 

The analysis described above was also applied 

to fractal structu res: the von Koch snowflakes 

(|Peitgen and JiirgenL 19921) and fractio nal Browian 
motion (fBm) images ( Turner et all Il998() . For the first 
object, whose fractal dimension is ~ 1.26, we produced 
snowflakes up to level 6 of different sizes (see Appendix) 
and studied their perimeter and area scaling. For fBms 
we created two sets of 150 images of expected fractal 
dimensions 1.8 and 1.6 respectively. Each f Bm image was 
segm ented with seven different thresholds llTurner et all 
Il998l) and perimeter and area of the structures selected 
by the different thresholds were combined to study the 
fractal dimension. 

In figEl the measured dimension dl is plotted as a 
function of object size for von Koch snowflakes of level 
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Fig. 6. dl evaluated for von Koch snowflakes of level 6. 
Like non fractal objects and real data, dl increases with 
object size and reaches a plateau at areas > 1000 pixeP. 
The plateau value, about 1.34, is an overestimate of the 
snowflake fractal dimension (see text). 

6. The fractal dimension increases with object size from 
a minimum value, in this case 1.15, to an almost con- 
stant value approximating the theoretical one, over the 
size range 1000 to 5000 pixel'^. The plot is surprisingly 
reminiscent of that found for real solar structures (§3). 
Note that the plateau value, about 1.34, exceeds the one 
expected theoretically. This reflects the overestimation of 
the snowflake perimeter inherent in the perimeter measure 
algorithm employed, as discussed previously for simple 
non fractal triangle. A rise of fractal dimension with ob- 
ject size was also observed for fBm images. Measurements 
of D are similarly affected by pixelization at small scales, 
with more complex structures harder to resolve and thus 
showing greater measurement error in the deduced fractal 
dimension, as shown in table El 

Both regular and fractal objects show similar pixeliza- 
tion induced errors in the fractal dimension estimation. 
These effects are greater at smaller areas, where the lack 
of resolution causes the objects to appear round, thus both 
dl and D increase rapidly with object size for object ar- 
eas less than ~ 500 — 1000 pixel^ and some hundred pixels 
square respectively. For objects of larger area, D and dl in- 
crease more slowly, but show deviations from the expected 
theoretical value reflecting how the structures map onto 
the pixel grid. 

We thus suggest that the minimum object size thresh- 
olds (about some tens of pixeP ) applied in previous works 
fe.g. lVogelaar fc Wakkeill994l) are insufficiently conserva- 
tive, with residual pixelization effects significantly infiu- 
encing the final results even for objects of ^ 1000 pixel^. 

5.2. Resolution and seeing effects 

The fractal dimension estimate of solar features depends 
on the resolution of the images analyzed. Resolution is de- 
termined not only by the detector pixel size (image scale), 
but also by the aperture of the telescope, any instrumen- 
tal aberration, and, for ground based instrumentation, the 
distortion introduced by atmospheric turbulence (seeing). 
Thus the pixel scale and the resolution are not the same. 
To evaluate the effects of resolution on the estimation of 



Dth D Dsm £ 

fBm 1.80 1.516 1.366 9.9% 
fBm 1.60 1.413 1.303 7.8% 
von Koch 1.26 1.310 1.273 2.8% 

Table 3. Fractal dimension measure for different frac- 
tals, before D and after Dsm smoothing by convolu- 
tion with a Gaussian of fwhm=2, and the relative error 
e = {D — Dsm)/D. Note that to distinguish resolution 
from pixelization induced effects, the error is evaluated 
respect to D and not to Dth- 
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Fig. 7. Facular fractal dimension estimated on the two 
different full resolution MLSO quality sets described in 
the text. When the estimation is carried out on images 
less affected by seeing degradation, the measured fractal 
dimension is higher. 



fractal dimension, we analyzed the scaling of dl and D 
with area after convolving von Koch snowflake and fBm 
images with Gaussian functions of different widths. We 
obtained, as one might have expected, a decrease in both 
dl and D accompanying the smoothing. Table Oof D (fit 
over the entire perimeter area range) shows also that the 
smoothing effects become more important as the structure 
complexity increases. 

A Gaussian function is a rough approximation to the 
seeing and instrumental aberration Point Spread Function 
in real images. Moreover, seeing is a time dependent phe- 
nomenon, so that images acquired at different times are 
affected by different degradation. In order to investigate 
directly the effect of variable seeing on the computed frac- 
tal dimension of structures in real data, we examined full 
resolution PSPT images from MLSO after selection based 
on the quality criteria described in §2.2. Images were seg- 
mented with the first technique explained in §2.3. FigureH 
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shows that a real decrease in resolution resulting from de- 
graded observing conditions leads to an underestimation 
of features' complexity at all scales. 

6. Interpretation of results measured 

We have shown that the determination of fractal dimen- 
sion of features in digitalized images by the two estima- 
tors dl and D is affected by pixelization and resolution. 
Understanding these effects is essential to the interpreta- 
tion of the results obtained from OAR and MLSO PSPT 
images (§3) as well as those previously reported from 
studies carried out with similar techniques on other data. 
Pixelization errors occur at all scales, but are generally 
more important for the smallest objects. This causes the 
estimated fractal dimension to increase rapidly with ob- 
ject size and become almost constant at areas larger than 
a critical threshold, that we estimated to be ~ 500 — 1000 
pixel^ for dl and some hundreds of pixeP for D. Seeing and 
instrumental induced image degradation smooths edges 
making structures appear rounder, resulting in a reduced 
fractal dimension. This effect is expected, on the basis 
of synthetic fractal data, to be more important for more 
complex objects (§5.2). 

We thus suggest that the rise of dl with object size 
observed in PSPT data, as well as for example in Meudon 
spectroheliograms esme- Ribes et al.l Il996|) and MDI 
magnetograms (|Meunier.. .1999^ i2004j) , is most likely an ef- 
fect of image pixelization, rather than a signature of an in- 
tr insic mu ltifractality of active regions. Conclusions drawn 
in lMeuni'en (^2004) concerning a change in physical proper- 
ties of magnetic structures at supergranular scales, should 
thus be reviewed in light of the results shown in this pa- 
per. Pixelization is als o likely the cause o f the 'break of 
similarity' observed bv l.Tanfien et al.l 1 20fl.ltl . the break oc- 
curring at the same pixel scale for both real and simulated 
images, in spite of the different physical scale implied, and 
in the same area range suggested by our synthetic fractal 
studies. 

For larger objects, the fractal dimension estimate is 
most affected by seeing. The values measured for large 
scale structures in the PSPT observations (1.5 to 1.7) 
are therefore likely an underestimate of the real value. 
Nevertheless, in some cases pixelization can cause over- 
estimation of the fractal dimension at largest areas, as 
shown for instance for the von Koch snowflake. We cannot 
therefore in principle exclude some compensation due to 
the combined effect of pixelization and reduction of res- 
olution. We finally note that MDI magnetograms, while 
not affected by seeing, are slightly defocu sed, so that the 
reso lution is twice that of the pixel scale JSch errer et al.L 
Il995i) . The same considerations made for results obtained 
with ground based measurements thus also apply to re- 
sults obtained with MDI. 

A study of the fractal geometry of solar active re- 
gions and its variation with the magnetic activity cycle 
is thus feasible if it focuses only on large features, em- 
ploys a constant segmentation technique throughout, and 



utilizes data of consistently high quality. The OAR PSPT 
images analyzed over a period of six years marginally meet 
these requirements. They do not show, however, a clear 
correlation with the solar cycle. Moreover the variations 
measured in D appear to be dominated b y variations in 
size d istribution of the examined features ijErmoUi et all 
which in turn weight the perimeter- area fit. 

7. Conclusions 

We have analyzed the fractal dimension of bright features 
identified in solar chromospheric Call K images obtained 
during the last six years, corresponding to the descending 
phase of solar cycle 23. The results obtained are in general 
agreement with those reported in literature, in particular 
with those studies that have been carried out with similar 
fractal estimators {dl and D). 

We have also investigated the effects of pixelization 
and resolution on the fractal dimension estimates ob- 
tained, studying these effects on real and simulated data, 
the latter including both non-fractal objects and ob- 
jects whose fractal properties are well known, von Koch 
snowflakes and fBm images. We have shown that fractal 
dimension estimates suffer from pixelization errors at all 
scales, but errors are generally more important for areas 
less than ^ 500 — 1000 pixeP. Particularly, pixelization 
causes measured fractal dimensions to increase with ob- 
ject size in the case of both fractal and non- fractal objects. 
Our results thus indicate that the increase of fractal di- 
mension with the feature size reported in literature by 
some previous analyses is likely an effect of pixelization 
and image degradation, rather than a signature of an in- 
trinsic multifractality of active regions. To reduce these 
effects, we thus suggest a restriction of the analyses to 
objects whose areas are larger then the quoted value. 

Our analyses also showed that image degradation due 
to both seeing and instrumental effects smooths features 
edges making them appear rounder. The fractal dimen- 
sion estimated is consequently lower than expected, even 
for large objects. Perhaps image de gradation effects ca. n 
be compensated for, as suggested in lJanBen et al.1 1)2003|) . 
Alternatively, a more careful estimate of these effects can 
be carried out by the analyses of images of known com- 
plexity objects (eg. fBms) convolved with realistic Point 
Spread Functions. We did not apply either of these com- 
pensations to our data and leave them for future investi- 
gation. 

Finally we note that, as demonstrated bvl Baveve et alJ 
(1998) for the box counting method, pixelization effects 
can influence other fractal dimension estimators as well, 
and careful quantitative measure of the effects for each 
measure employed is essential to the interpretation of the 
results. We have not addressed this problem in this work 
and leave also this issue for future research. 
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Appendix A: On the estimation of D and dl on 
digitalized images 

The estimators adopted in this paper to evaluate frac- 
tal dimension of selected structures are based on the 
perimeter-area relation. This consists of measuring the 
perimeter and area of a structure at different resolutions. 
For regular structures, perimeter scales as the square root 
of the area, while for fractal structures the exponent is 
greater than 0.5. By definition, the exponent is the fractal 
dimension of the studied object. Note that in §2.5 we nor- 
malized the exponent so that the fractal dimension is one 
for regular structures and greater than one for fractal ob- 
jects. It's worth noticing that, when adopting this method 
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Fig. A.l. Perimeter (in units of pixel) and Area (in 
units of pixels square) in logarithmic scale of detected 
structures on OAR PSPT data taken during summer 
2002. Continuous line is the fit to the whole set of data 
(_D=1.35). Points at area larger than about 1000 pixels 
square are better approximated by a higher slope line. 
Horizontal line is the area window width over which dl is 
estimated. 



to investigate fractal nature of magnetic solar regions, one 
implicitly assumes that the different size selected struc- 
tures are the 'same' object observed at different resolu- 
tion. 

When plotting the perimeter-area relation of real data, 
one expects to find at least three different regimes. At the 
smallest sizes, because of resolution, an object's detail is is 
not fully detected, thus perimeter and area scale as for reg- 
ular non-fractal structures. At the largest areas, a break 
in similarity can occur for physical reasons (for instance 
the object under study has a finite maximum size above 
which it is no longer a fractal). At intermediate areas the 
object scales as a fractal. As an example, in fig. lA.ll we 
plot the perimeter-area relation for structures selected in 
OAR-PSPT images. The straight line is a linear fit to the 
points. It's evident that points don't lie on a single line, 
but rather on a curve. One is thus tempted to measure 
the tangent to the curve as a 'local' measure of the fractal 
dimension. This is what the measure in dl employed in 
the text does. All our plots (cfr. figures [?I3I51 1617(1 showed 
that dl increases with object size, becoming eventually 
constant at some typical scale (generally in the area range 
500-1000 square pixels). The measured dl, as well as the 
area at which it becomes constant, are functions of the 
window size. Nevertheless, when using this estimator, it 
is important to keep in mind that any fractal estimate re- 
quires the autosimilarity to be valid over some orders of 
magnitude ijBaveve et al.L Il998|) (the area range A log A 
=1.5 adopted for dl in this and in other works is therefore 
in principle too small). 
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Fig. A. 2. Left: Perimeter versus area in logarithmic scale of snowflakes of levels 2,4, and 6. Because of pixelization, 
these structures exhibit a fractal scaling only within certain area ranges, the bounds depending on the snowflake level. 
Right: dl versus area evaluated with different window sizes for snowflakes of level 6. Peaks obtained with the small 
window (open circles) are due to the steep variations visible in plot on the left. Peaks are not detected with a larger 
window (full circles). The area range over which dl is almost constant is the range over which the simulated images 
are fractal. 



The results we obtained with von Koch snowflakes il- 
lustrate clearly these issues, von Koch snowflake images 
of differe nt sizes were prodiiced fo llowing the iterative 
scheme of lPeitgen and JiireeiJ l)l992() . After each iteration, 
or level, the snowflake is more structured, with an increase 
in both perimeter and area. In the limit of infinite itera- 
tions, the perimeter tends toward infinity and the area 
approaches a finite value. Here we investigate structures 
constructed with up to 6 levels. The fractal dimension of 
the von Koch snowflake is log 4/ log 3 ^ 1.26. 
In left panel of fig. I A. 21 the perimeter-area relationships 
for snowflakes of levels 2, 4, and 6 are plotted with log- 
arithmic scaling. For each level, the relationship traces a 
curve made up of segments whose slope is 1/2 connected 
by segments of slope greater than 1/2. At largest areas 
all the points lay on parallel lines of slope 1/2. At those 
scales the snowflakes of all the represented levels are fully 
resolved on the grid employed. As the dimensions of the 
objects are reduced, fewer details at any fixed construction 
level are resolved, the measured perimeter decreases at a 
rate faster than A^^'^, and the perimeter-area curve steep- 
ens. The slope flattens to a value of one-half again each 
time the grid resolution is sufficient to capture the details 
of the next lower level. Finally, at smallest areas most ge- 
ometric details are lost and all the objects, independent 
of their initial construction level, appear non-fractal. 

The scaling of dl better reflects the change in slope 
with objects size. As an example, in rig ht panel of figE21 
we show results obtained for level 6. Here full and open 
dots represent respectively dl obtained with a window 
of Alog^ =1.5 and a window of Alog^ = 0.5. With 
the largest window only the slope change that occurs at 
largest areas is visible. The others occur on scales smaller 
than the window so that they are not 'detected' and a 
plateau is observed. At smallest areas dl drops because 
of the resolution effects explained before. When a smaller 
window is used, 6 peaks are visible, corresponding to the 
6 slope-changes visible in the perimeter area scatter plot. 
In this case, there is an area range over which fractal 
dimension oscillates around a constant value. At smallest 



areas larger amplitude oscillations are observed. Both 
curves show clearly two of the three regimes mentioned 
above. The object scales as a fractal in the range 3 
< log(Area) < 4.5. At smaller areas, pixelization effects 
dominate the measurements because resolution is insuffi- 
cient to allow detection of all the structure's details. The 
third regime is not evident in the two curves, since not 
enough points are available at largest areas to perform 
the fit. If more points were available we would observe a 
decrease of dl toward the value of one, as slightly visible 
at largest areas for fits performed on the smaller window. 



